Extended wetlands priority habitats

from pathlib import Path
import pandas as pd

import os
os.environ['USE_PYGEOS'] = '0'
import geopandas as gpd


# import colorcet as cc
import datashader as ds
# import datashader.transfer_functions as tf

import holoviews as hv
from holoviews.operation.datashader import rasterize
from holoviews import opts

import geoviews as gv
import cartopy.crs as ccrs

from spatialpandas import GeoDataFrame

from dask import dataframe as dd

gv.extension('bokeh')

• Prepared plots showing extended wetland priority habitats.

• Tried Holoviews/Geoviews/Datashader for interactive plots, they works very well.

• Tried jupyterlite and pyodide but there are no wheels available for datashader and co. so micropip was not able to run with them

import requests
import io

import sys
if not sys.warnoptions:
    import warnings
    warnings.simplefilter("ignore")

url = 'https://github.com/iacopoff/iacopoff.github.io/raw/main/data/n2k_extwet_priority_poly_onehot_s.fgb'

dn = requests.get(url).content

df = gpd.read_file(io.BytesIO(dn))

url2 = 'https://github.com/iacopoff/iacopoff.github.io/raw/main/data/annex_I_extended_wetlands.csv'
dn2 = requests.get(url2).content
lookup = pd.read_csv(io.StringIO(dn2.decode('utf-8'))).set_index('habitatcode')

# dpath = Path('~').expanduser()  / 'eos_data/KCEO/uc_wetlands/n2k_aoi/GIS'

#fpath = dpath / 'n2k_extwet_priority_poly_onehot_s.fgb'

#df = gpd.read_file(fpath)
#lookup = pd.read_csv(dpath / 'annex_I_extended_wetlands.csv').set_index('habitatcode')

hcodes = df.filter(regex='\d',axis=1).columns.tolist()

df['n_pr_hb'] = df.filter(regex='\d',axis=1).sum(axis=1)

df = df.to_crs(3857)

df_sp = GeoDataFrame(df)

ddf = dd.from_pandas(df_sp, npartitions=10)

# cvs = ds.Canvas()
# tf.shade(cvs.polygons(ddf[ddf['3170'] == 1], geometry='geometry', agg=ds.mean('n_pr_hb')), cmap=cc.kg);

bm = hv.element.tiles.EsriUSATopo()

def get_canvas(hcodes, bm):
    out = []
    for h in hcodes:
        out.append(bm * rasterize(hv.Polygons(ddf[ddf[h] == 1], vdims=['sitename','n_pr_hb']), aggregator=ds.mean('n_pr_hb')).opts(
            colorbar=True, cmap='bgy',fontsize={'title':8}))
    return out

outs = get_canvas(hcodes, bm)

# def compute_partitions(el):
#     n = ddf[ddf['3170'] == 1].cx_partitions[slice(*el.range('x')), slice(*el.range('y'))].npartitions
#     return el.opts(title=f'Population by country (npartitions: {n})')

# bm * out1.apply(compute_partitions).opts(width=700, height=500, tools=["hover"])

# opts.Polygons(logz=True, tools=['hover'], xaxis=None, yaxis=None,
#                show_grid=False, show_frame=False, width=500, height=500,
#                color_index='das', colorbar=True, toolbar='above', line_color='white')

l = hv.Layout([ts.opts(width=400, height=300).relabel(name + ": " + lookup.loc[name][0]) for name, ts in zip(hcodes,outs)]).opts(
    opts.Tiles(xaxis=None, yaxis=None, width=100, height=100)).cols(3)

l

# import matplotlib.pyplot as plt
# import cartopy.crs as ccrs

# from cartopy.io.img_tiles import Stamen


# tiler = Stamen('terrain-background')
# mercator = tiler.crs

# # df = df.set_index('sitecode')

# def plot_fig(hcodes):
#     for h in hcodes:
#         fig = plt.figure(figsize=(10,8))
#         ax = fig.add_subplot(1, 1, 1, projection=mercator)
#         ax.set_extent([-15, 37, 34, 71], crs=ccrs.PlateCarree())

#         ax.add_image(tiler, 6, alpha=0.8)
#         ax.coastlines('10m',color='grey', )
#         df[df[h] ==1].plot(ax=ax, color='red', linewidth=0.5) # column='n_pr_hb', ax=ax, categorical=True, legend=True, cmap='red')

        
#         plt.title(h + ': ' + lookup.loc[h][0])
#         fig.savefig(dpath / f'../habitats/habitat_{h}.png', dpi=300)

#plot_fig(hcodes)

# fig = plt.figure(figsize=(15,13))
# ax = fig.add_subplot(1, 1, 1, projection=mercator)
# ax.set_extent([-15, 38, 34, 72], crs=ccrs.PlateCarree())

# ax.add_image(tiler, 6, alpha=0.8)

# df.plot(column='n_pr_hb', ax=ax, categorical=True, legend=True, cmap='jet')

# #ax.coastlines('10m',color='grey')
# plt.title("Count of 'wet' priority habitats in N2K sites")
# fig.savefig(dpath / f'../habitats/priority_habitats_count.png', dpi=300)